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Resolving the CO Snow Line in the Disk around HD 163296 
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ABSTRACT 



We report Submillimeter Array (SMA) observations of CO (J=2-l, 3-2 and 
6-5) and its isotopologues (^=^C0 J=2-l, C^^O J=2-l and C^^O J=3-2) in the 
disk around the Herbig Ae star HD 163296 at ~ 2" (250 AU) resolution, and 
interpret these data in the framework of a model that constrains the radial and 
P^ ' vertical location of the line emission regions. First, we develop a physically self- 

^ ■ consistent accretion disk model with an exponentially tapered edge that matches 

"-^ ■ the spectral energy distribution and spatially resolved millimeter dust continuum 

I ' emission. Then, we refine the vertical structure of the model using wide range of 

i-{ • excitation conditions sampled by the CO lines, in particular the rarely observed 

52 ■ J=6-5 transition. By fitting ^^CO data in this structure, we further constrain 

the vertical distribution of CO to lie between a lower boundary below which CO 
freezes out onto dust grains (T < 19 K) and an upper boundary above which 
CO can be photodissociated (the hydrogen column density from the disk surface 
^ ■ is < 10^^ cm~^). The freeze-out at 19 K leads to a significant drop in the gas- 

i/-^ ■ phase CO column density beyond a radius of ~155 AU, a "CO snow line" that we 

[■^ ■ directly resolve. By fitting the abundances of all CO isotopologues, we derive iso- 

^-^ ■ topic ratios of ^^C/^^C, ^^O/^^O and ^^O/^'^'O that are consistent with quiescent 

interstellar gas-phase values. This detailed model of the HD 163296 disk demon- 
strates the potential of a staged, parametric technique for constructing unified 
^ . gas and dust structure models and constraining the distribution of molecular 

j^ ■ abundances using resolved multi-transition, multi-isotope observations. 
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Introduction 



The disks around pre-main sequence stars are the reservoirs of raw material that repre- 
sent the initial conditions for the formation of planetary systems. The spatial distribution 
of mass in these disks is a fundamental property, as it sets a constraint on the contents and 
orbital architectures of planetary systems. In that sense, measurements of disk densities and 
temperatures can provide strong, albeit complex, constraints on planet-forming scenarios. 
Because the dominant constituent of these disks is believed to be cold H2 gas that is gen- 
erally unobservable, our knowledge of many disk properties (including mass) relies on the 
observation and interpretation of minor constituents. Dust has naturally been the focus of 
the most work on disk str ucture and evolution, as it dominates the opacity and emits bright 
continuum radiation (e.g.. [Andrews et al.ll2009l . I2OIOI ). However, converting those continuum 



surface brightnesses to disk (gas) densities is challenging due to the large uncertainties in 
the assumed dust properties (composition, size distribution, etc.) and the dust abundance 
relative to the gas, which may deviate substantially from the canonical interstellar value 
(0.01) and vary spatially in the disk. In principle, spectral line emission from trace gas 
species can provide complementary constraints on the disk structure, but the interpretation 
of the observations can also be difficult. These emission lines depend on the gas abundances 
that are set by chemical and physical processes, as well as the excitation conditions set by 
the local densities, temperatures, and incident radiation fields. With a dataset of sufficient 
quality, these complexiti es can be leverage d into powerful probes of otherwise inaccessible 



disk characteristics (e.g.. lKamp et al.ll2010l ). 



The chemical property of CO, the most abundant molecule after H2, is thought to be 
well understood in disks. The spatial distribution of gas-phase CO is expected to be shaped 
by photodissociation both near the star and high in the disk atmosphere, as well as by 
a depletion process when it freezes onto dust grair is at low ternperatures (below ^20 K) 
at large radii and deep in the disk interior (e.g., lAikawa et al.l Il996t lAikawa fc Nomura 
20061 : iGorti fc Hollenbachll2008l ). The peak abundance of CO is also pre dicted to be almost 
independent of radius, around 10~^, the typical value in molecular clouds (lAikawa fc Nomura 
20061 ). However, observational studies of CO emission from disks have yielded puzzling 
results. In general, the CO in disks is found to be under-abundant by factors of 10-100 
compared to molecular clouds, a deficit usually attributed to CO depletion in cold, dense gas 
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flDutrev et al.lll994 Il996l . Il997l : Ivan Zadelhoff et al.ll200ll : iThi et al.ll2004h . However, multi- 
transition and multi-isotope CO studies report a wid e range of relative depletions of CO and 
its isotopologues. For example, iDartois et al.l ( 120031 ) fin ds that all o f the C O isotopologues 
are depleted by an identical uniform factor of ~10, but iPietu et al.l ( 120071 ) finds in a small 
sample of disks that the ^^CO/^^CO ratios in the outer part of disks are much lower than 
the standard ^^C/^^C ratio in the solar neighborhood: these low values were attributed to 
significant carbon fractionation. At the same time, both of these studies suggested that the 
outer radius derived from ^'^CO is smaller than the one derived from ^^CO indicating their 
^^CO/^'^CO ratios are abnormally large at the outer disk edge of their models, which leads 
to speculations on the imp ortance of selective photo-dissociation due to self-shielding in the 
photolysis of gaseous CO (JDutrey et al.l 120071 ). They also find that a large amount of CO 
remains gaseous at temperatures as low as 10 K. Possible e xplanations for the presence of 
the inf erred cold CO inclu de vertical mixing (JAikawal 120071 ) and photodesorption off dust 
grains (JHersant et al.l 120091 ). All above depict a very complex picture of the distribution of 
CO gas in the disks, which remains to be explained by detailed chemical models combined 
with realistic physical disk structures. However, it is unclear, how much of these variations 
depend on the individual sources, and how much, if any, depends on the modeling approach. 

A robust methodology for interpreting CO observations is a high priority for our efforts 
to constrain the total gas distribution in disks, and to provide a basis for interpreting ob- 
servations of more complex molecules. In this study, we develop a modeling framework that 
is able to account for observations of both the dust continuum and CO line emission from a 
protoplanetary disk. Specifically, we combine a previously-established dust emission model- 
ing formalism with multiple transitions, spatially resolved CO line data (CO J=2-l, 3-2 and 
6-5, i^CO J=2-l, C^^O J=2-l and C^^O J=3-2) to construct a simplified and internally 
self-consi stent model of the dis k around the Herbig Ae star HD 163296. At a distance of 
~122 pc (jPerryman et al.lll997l ). this 2.3 M© star (spectral type Al; age ~4 Myr) harbors a 



large disk with strong m illimeter continuum and molecular line emission (IOill200ll : iThi et al. 



20041 : iNatta et al.ll2004l ). It does not seem to be associated with any known star- formation 
region, dark clouds, or reflection nebulae. The disk around HD 163296 exhibits a scat- 



tered light pattern extending out to a radius of ~500 AU (JGrady et al.ll2000l ). perpendicular 
to a bipolar microj e t traced by a chain of HH knots visib le in coronagraphic observations 
(JDevine et al.l I2OOOI : ICradv et al.l I2OOOI : IWassell et al.l I2OO6I ) . Millimeter and submillimeter 
interferometric observations have explored this large disk an d its Keplerian velocity field in 
low- J CO lines (JMannings fc Sargentlll997l : llsella et al.ll2007l ). The large disk size and strong 
molecular emission make this system an excellent target for resolved millimeter observations 
and detailed modeling. 



This paper is organized as follows. We describe the CO and mm continuum data in §2, 
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and present the main observational results in §3. Our extensive modeling effort is detailed 
in §4. We discuss the modeling results and compare with previous work on this subject in 
§5, and provide a brief summary in §6. 



2. Observations 

Observations of HD 163296 (RA: 17^^56"21!279, DEC: -21°57'22'.'09; J2000.0) were con- 
ducted between 2005 August and 2010 September using the 8-antenna SMAl^ interferometer 
located atop Mauna Kea, Hawaii. Table 1 provides a general summary of the observational 
parameters. 

For the ob servations of the CO 6-5 line, we used the correlator settings adopted by 



Oi et al.l (120061 ) in their similar study of the TW Hya disk. At 690 GHz, there is no nearby 
quasar bright enough to use for phase referencing with the SMA. However, Callisto was 
located only 10° away from HD 163296 during our observation on 2007 March 20, and was 
monitored every 20 minutes for use in the gain and absolute flux calibration. At that time, 
Callisto had a diameter of 1'.'3 and a zero-spacing flux density of 48.1 Jy at 690 GHz. Based 
on the uncertainties of the Callisto emission model, we estimate a 10% systematic uncertainty 
in this adopted flux scale. The 690 GHz continuum emission from HD 163296 was sufficiently 
bright (with a flux density of 7.5 Jy) that we could correct the gain response of the SMA 
with a single phase-only self-calibration iteration. 

The observations of the CO 2-1 line took place on 2010 May 17 in the compact conflgu- 
ration and 2010 September 14 in the extended conflguration. The correlator was conflgured 
to include 2048 channels in the 104 MHz segment of the correlator centered on the rest 
frequency of the line. The remaining 1.3 GHz of the correlator bandwidth in both side- 
bands were conflgured with a uniform spectral resolution of 256 channels in each 104 MHz 
correlator segment to achieve maximum continuum sensitivity. The weather was good on 
both nights, with the 225 GHz opacity just below 0.1 and stable atmospheric phase. The 
nearby quasar 1733-130 (10.5 degrees from the source) was used to cahbrate the atmospheric 
and instrumental phase, and the solution was checked using the quasar 1911-201, which was 
included for a short period during each loop between source and calibrator. Uranus and 
Neptune were used as flux calibrators on September 14, yielding a 1733-130 flux of 2.71 Jy; 
Callisto was used to calibrate the May 17 flux, yielding a 1733-130 flux of 2.04 Jy. 



^The Submillimeter Array is a joint project between the Smithsonian Astrophysical Observatory and the 
Academia Sinica Institute of Astronomy and Astrophysics, and is funded by the Smithsonian Institution and 
the Academia Sinica. 
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The observations of ^^CO 2-1 and C^^O 2-1 were carried out simultaneously on 2010 
May 15 in the compact configuration and 2010 September 11 in the extended configuration. 
The correlator configuration was identical to the CO 2-1 setup, except that the high reso- 
lution correlator chunk was centered on the rest frequency of ^^CO 2-1 with C^^O 2-1 in 
a nearby low-resolution chunk. The remaining 1.2 GHz were devoted to uniform spectral 
resolution continuum observations. The weather was excellent on both nights with stable 
phase and 225 GHz opacity varying between 0.04 (primarily during the extended track) and 
0.06. The same calibrators as for the C0(2-l) observations were used; the derived flux of 
1733-130 was 2.35 on September 11 and 2.04 Jy on May 15. 

The observation of the C^^O 3-2 line took place on 2005 August 23 in the compact con- 
figuration. The correlator was configured with a uniform spectral resolution of 128 channels 
over 104 MHz, which provided 0.8 MHz frequency resolution. The nearby quasars 1833-210 
and 1921-293 were used as the calibrators. Uranus was used as flux calibrator and the derived 
fluxes of 1833-210 and 1921-293 were 0.57 and 3.76 Jy, respectively. 



The observations of the CO 3-2 observations were summarized i n l Hughes et al.l ( 1201 ll ). 
All of the calibration was performed using the MIR software package!. Images of the con- 
tinuum and the spectral lines were generated and CLEANed using standard techniques in 
the MIRIAD software package. 



3. Results 

Figure [U shows the disk-averaged spectra of the CO J=2-l, 3-2 and 6-5 lines, along 
with the J=2-l lines of ^^CO and C^^O and the J=3-2 line of C^^O 3-2. These spectra were 
extracted from the SMA channel maps in 10" square boxes centered on HD 163296, except 
for the CO 6-5 line. That fainter spectrum was produced from a 4"-wide box chosen to cover 
the area with statistically significant emission. Each emission line shows the double-peaked 
velocity profile characteristic of disk rotation. The low-J CO lines have an additional feature 
in their spectra at Vlsr = 13 km s~^, which is distinct ly offset from the disk kinematically 



and much more obvious in previous single-dish spectra (JThi et al.ll200ll . |2004J ). That feature 
is apparently spatially extended, and likely is associated with foreground and/or background 
molecular clouds. Table 2 lists the integrated intensities of the observed lines and the con- 
tinuum flux densities. Note that the signal-to-noise ratios range from over 200 for the CO 
2-1 line to ~10 for the CO 6-5 and C^^O 3-2 lines. Figure 121 [3] an d HI show the c hanne l 



maps of all the transitions except for CO 3-2 which were presented in lHughes et al.l ( 120111 ). 



^http://www.cfa. harvard. edu/'^cqi/mircook. html 
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4. Analysis 

As a first step in developing a framework for understanding the cfiemical structures of 
protoplanetary disks, we aim to interpret this suite of observations in the context of an inter- 
nally self-consistent model for the physical conditions of the gas and dust in the HD 163296 
disk. Our modeling analysis consists of three distinct stages. First, we construct a tem- 
plate two-dimensional model for the disk structure based on observations of the dust, using 
the broadband spectral energy distribution (SED) and the resolved millimeter continuum 
emission as diagnostics (§4.1). Second, we refine the vertical temperature distribution of 
that template structure based on the spatially resolved line ratios of the optically thick CO 
transitions (§4.2). And third, we constrain the radial and vertical CO abundance pattern 
using the emission from the more optically thin CO isotopologue transitions (§4.3). 



4.1. Stage 1: Template Disk Structure Model 



Our starting point is with a set of disk structures calculated following the presc ription for 
a steady viscous accretion disk model (see iDAlessio et al.lll998l Il999l . 1200 ll . 120061 ). In these 
models, the gas surface density is determined based on the conservatior i of angular momen- 
tum flux and depends on the mass flow rate (M), viscosity coefficient a (IShakura fc Sunyaev 
19731 ). and midplane temperature (Tq), such that S oc M jaT^. We assume that M is constant 
throughout the disk, but make an important structural modification from previous genera- 
tions of these models. Motiva ted by the similarity solution for the time e volution of accretio n 
disks (JHartmann et al.lll998[ ) and recent millimeter- wave observations ( iHughes et al.ll2008l ). 
we allow the viscosity coefficient to vary radially such that a = a^ exp R/Rc. This change 
makes no difference at small radii {R <^ Re), but effectively adds an exponential taper to 
the surface density profile outside the characteristic radius. Re- Although this modification 
mimics the E profile derived from the similarity solution models, it does not reproduce the 
predicted behavior of the mass fiow rate as a function of radius. However, that M profile is 
not expected to have any observable effects on the disk structure at these large radii, since 
stellar irradiation is a much more important heating mechanism than viscous dissipation in 
those regions. 

For a given flow rate (M), viscosity coefficient (ao), and characteristic radius (Rr), the 



density and temperature structure of this model is determined as described by lD'Alessio et al. 



( 1l998l . 119991 ) . We consider heating from the mechanical work of viscous dissipation (relevant 
only in the midplane of the inner disk), accretion shocks at the stellar surface, and passive 
stellar irradiation, and follow the radiative transfer of that energy with I+ID calculations 
using the Eddington approximation and a set of mean dust opacities (gas opacities are con- 
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sidered negligible). The dust is assumed to be a mixture of segregated spheres composed of 
"astronomical" silicates and graphite, with abun dances (relative to the total gas mass) of 
Csii = 0.004 and Cgra = 0.0025 (JDraine fc Led 1 19841 ): the "reference" dust-to-gas mass ratio is 
Cref = 0.0065. At any given location in the disk, the grain size (a) distribution of these dust 
particles is assumed to be a power-law, n{a) oc a~^-^, between amin = 0.005 /xm and a spec- 
ified Omax- We assume the disk has two grain populations, each with a different maximum 
size. The "small" grains are distributed in the disk atmosphere and have amax = 0.25 /xm 
(as in the interstellar medium), and the "big" grains are concen trated toward the midp lane 
and have Omax = 1 mm. Utilizing the dust settling prescription of iD'Alessio et al.l ( 120061 ) . the 
"small" grain population in the upper disk layers is assumed to have a dust-to-gas mass ratio 
eCrcf and the "big" grain population near the midplane has a dust-to-gas ratio determined so 
that mass is conserved at each radius. A smooth vertical transition between the two grain 
populations is made at a height zug above the midplane. 

The inner boundary of the structure model is taken to be the radius where silicates are 
sublimated (whe re T ^ 1500 K). B e cause this inner rim is irradiated by the photosphere of 
the central sta r (INatta et al.ll2000l: iDuUemond et al.ll200ll ) and the accretion shock at the 



stellar surface ( MuzeroUe et al.l 120041 ) at a normal incidence angle, the material is heated to 



high temperatures and therefore extends (or "puffs" up) to a larger vertical height, Zwaii- The 
ra dial structure of that "w all" feature (and its emergent spectrum) is calculated as described 
by ID'Alessio et al.l (|2005l ). The structure model is truncated at 540 AU, large enough to 
have no effect on our analysis. 

In the model cal c ulatio ns, we adopt the HD 163296 stellar parameters advocated by 
van den Ancker et al.l (119981 ) - M* = 2.3 M©, R^, = 2Rq, and T* = 9333 K - and a flow rate 
equiv alent to the accretion rate onto the s tar, M = 7.6 x 10~ ^ M© yr~^ (JGarcia Lopez et al. 



20061 ). Based on the CO observations of llsella et al.l (120071 ). we fixed the disk inclination 
angle to i = 44° and the major axis position angle to 133° east of north. For this set of fixed 
parameters, the inner wall is located at a radius of 0.6 AU. There are five remaining free 
parameters in the model: {ao, Re, e, 2;waii, ^big}- To identify a model structure that provides 
a reasonable match to the HD 163296 observations, we started by comparing synthetic data 
for a coarse grid of these parameters with the broadband SED. That parameter search was 
then refined by comparing the model predictions with the observed millimeter continuum 
visibilities at 271 and 341 GHz. Although not optimized for robust parameter estimation, this 
approach yielded a template model that exhibits a satisfactory match with these diagnostics 
of the dust disk. The adopted model has ao = 0.019, Re = 150 AU, e = 0.003 (the dust- 
to-gas ratio in the atmosphere is 0.3% of the reference value, Cref), and z„aii = 0.1 AU. The 
total mass of the model is 0.089 M©. All the stellar and disk properties are summarized 
in Table [31 Figure [5] shows how the value of the Re parameter affects the corresponding 



continuum visibility profiles at 271 GHz and 341 GHz, the frequencies with the best resolved 
millimeter continuum data. We find the visibility profi l es are matched best with R^ at 150 
AU, a bit larger than 125 AU derived by iHughes et al.l ( 120081 ) with the similarity solutions. 



Figure [6] compares the HD 163296 SED with the models with different values of ^big, 
ranging from 0.5H to 2.5/7, where H is the pressure scale height calculated from the local 
temperature value {H = c^/fi, where Cg is the sound speed and Q is the Keplerian angular 
velocity of the gas. Note that this is merely a convenient scaling, as the vertical structure 
is solved consistently and does not assume a vertically isothermal profile). The different 
behaviors of the models are likely caused by the change in the shape of the irradiation surface, 
which determines the fraction of stellar radiative flux intercepted and reprocessed by the disk. 
Apparently it is difficult to constrain Zbig from the SED alone due to the complexity of the 
SED modeling in the mid and far-IR wavelengths. The most substantial discrepancy between 
the observed SED and the model predictions can be noted in the strength of the 10 and 18 
/im silicate emission bands. This might be a reflection that the upper atmosphere layers are 
hotter than in the model, or that the optical properties of the grains in the atmosphere are 
slightly different (in size distribution and/or composition) than we assume. The modeling 
formalism we are using assumes the gas and dust are in thermal equilibrium at all locations 
in the disk, even in the upper layers of the atmosphere. While that is expected to be a good 
approximation at the height where most of the stellar radiation is deposited, and where most 
of th e emission from the bands is produced, the gas at larger heights may be even hotter 



e.g., Glassgold et al. 



2004 



Karnp &: DuUemond 



2004 



Jonkheid et al. 2004; Nomura et al. 



20071 : iGorti fc Hollenbachll2008l : IWoitke et al.l 120091 : iKamp et al.l 120101 ) . which might change 
the dust density distribution, and hence the silicate emission. 



4.2. Stage 2: Refining the Vertical Structure 

At the large disk radii probed by our CO observations, the vertical structure of the disk 
is determined by irradiation. The small dust grains in the disk atmosphere absorb energy 
from the incident stellar radiation field, and some of that energy is then re-emitted down 
into the disk interior at longer, infrared wavelengths. Given the increasing densities near the 
disk midplane, this "external" heating of the disk surface naturally produces a s tructure with 



a ver t ical temperature inversiori : the midplane is cool er than the atmosphere (ICalvet et al. 



199ll : IChiang fc Goldreichlll997l : iD'Alessio et al.lll998l ). In our model prescription for dust 
settling, the deeper layers in the disk are populated by big dust grains that have low infrared 
opacities, and therefore are heated less efficiently than their smaller counterparts. Therefore, 
this concentration of big grains near the disk midplane actually amplifies the temperature 
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contrast between the surface and interior. In practice, the location of the transition between 
the small and big grains, z^g, has a pronounced effect on the vertical temperature gradient 
(and therefore the vertical density structure). In Figure [71 we demonstrate how 2;big impacts 
the vertical distributions of temperatures and densities at a fixed radius of 200 AU. A more 
condensed population of big grains (lower zug) permits more heating at deeper depths into 
the disk (i.e., lower z), producing a warmer disk interior. Likewise, a more vertically extended 
population of big grains (higher zug) produces much lower temperatures in the disk interior. 

The effects of ^big are difficult to distinguish from the dust tracers alone: the infrared 
SED and millimeter visibilities do not effectively probe the shape of the vertical temperature 
profile. However, the CO line emission is expected to be generated in an intermediate disk 
layer that should be sensitive to this temperature inversion. The main isotope CO lines are 
optically thick, and therefore excellent temperature diagnostics. By measuring the emission 
from several CO rotational transitions, their relative strengths may be used to trace the 
temperature at different depths in the disk (each line probing a layer commensurate with 
its excitation). Previous analysis of multi-transiti on CO data have be en used to successfully 
measure the vertical temperature gradient (e.g., iDartois et al.ll2003l ). In this stage of the 
modeling, we utilized the resolved CO 2-1, 3-2, and 6-5 emission lines toward the HD 163296 
disk to determine the value of Zbig that best reproduces the inferred temperature gradient 
present in its intermediate layer. 

To compare the model predictions with the CO data, we first need to fix a set of geo- 
metric and kinematic parameters that affect the observed spatio-kinematic behavior of the 
disk. We assume the disk material orbits the central star in Keplerian motion, and fix 
the stellar ma ss and position, and non-thermal turbulent velocity width {dVturb) based on 
the models of iHughes et al.l (120111 ) and the references noted in §4.1. The effects of inclina- 
tion (i), position angle of the major axis (PA), and systemic velocity (Visr) are essentially 
orthogonal to those that depend on the assumed disk structure, and therefore can be op- 
timized independently. All the disk geometric and kinematic parameters are summarized 
in Table [31 We further assume a simple CO abundance model that dep ends on the local 
density and temperature via the parameters introduced by lOi et al.l (120081 ) . described below. 
For a given disk structure and CO abundance model, we use the non-local thermodynamic 
equilibrium two-dimensional accelerated Monte-Carlo radiative transfer code RATRAN to cal- 
culate the molecular excitation an d generate a sky-projected set of synthetic CO data cubes 
(JHogerheijde fc van der Takll2000l ). sampled at the same spatial frequencies as each SMA 
dataset. The coUisional rates are taken from the Leiden Atomic and Molecular Database 



(ISchoier et al.l |2005[ ) for non-LTE line radiative transfer calculations. Specifically, in the 
models used i n this pape r , we h ave used the new set of the CO collisional rate coefficients 
calculated by lYang et al.l (120101 ) . 
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Our adopted abundance model was introduced by iQi et al.l (120081 ) . and assumes that 
the CO emission originates in a vertical layer of the disk with a constant abundance. The 
upper (surface) and lower (midplane) boundaries of that layer are defined by the parameters 
as and 0"^, which represent vertically- integrated hydrogen column densities from the disk 
surface in units of 1.59 x 10^^ cm~^ (the conversion factor of the hydrogen column to A^ for 
interstellar dust). The CO abundance (/co) is assumed to be a constant for the hydrogen 
column densities from the disk surface cr^ > iV > cr^, with respect to the H2 density specified 
in the disk structure model. 

We compute model CO visibilities for a range of disk structure models with different 
2^big values, ranging from 0.5-3.0if in 0.5if steps, where H is the pressure scale height as 
we described before. For each Zbig value, we optimized the CO abundance model (o"m, cr^ 
and fractional abundance /co) to minimize the combined x^ value in reference to the CO 
J=2-l and 3-2 visibility data. Because the sensitivity of our CO J=6-5 data is relatively 
modest, this high-excitation line is better suited to an a posteriori check on the fitting results. 
Nevertheless, the J=6-5 line has a sufficiently high excitation that its luminosity relative to 
the low-lying J=2-l and 3-2 lines provides perhaps the best discriminant between different 
2;big values. Figure [S] shows the results of this stage of the modeling analysis, marking a 
direct comparison between the observed and model spectra of the main CO transitions for 
3 different ^big values - 1.5, 2.0, and 2.5if . While any of these models provides a suitable 
match to both the J=2-l and 3-2 spectra, only the Zbig = 2-ff model also matches the J=6-5 
line. At R=100 AU of this model, the abundance of big grains has decreased to 50% of its 
maximum value at a height around 16 AU from the midplane. Figure shows our favored 
two-dimensional density and temperature distributions for this 2;big = 2-ff model. 



4.3. Stage 3: CO Abundance Distributions 

We have established a gas-|-dust structure model that reproduces well the broadband 
SED, resolved millimeter continuum images, and a multi-transition set of CO spectral im- 
ages. The next step in the modeling analysis is to constrain the spatial distribution of the 
CO abundance, relying on the spatially resolved - and more optically thin - CO isotopologue 
emission, which probe much deeper into the midplane. Here we assume CO and its isotopo- 
logues share the same spatial distribution and only differ in fractional abundances. Instead 
of (Tm, at this stage we adopt the lower boundary (toward midplane) governed by the CO 
freeze-out temperature, Tco, such that /co -^ when T < Tco, which provides a way of 
interpreting this lower boundary more physically. For the upper boundary (toward surface), 
we still adopt ag as before. 
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Using the fixed structure model derived in §4.1 and 4.2, we compute a grid of syntlietic 
^^CO J=2-l visibility datasets over a range of a^, Tqq, and fisco values and compare with 
the observations. Figure [TD] shows the x^ surfaces for the ^^CO J=2-l emission in the space 
of these three parameters. We find that the data are best described when a^ = 0.79 ± 0.03 
(i.e., the hydrogen column densities from the disk surface are below 1.2-1.3x10^^ cm~^), 
Tco = 19.0 ± 0.3 K, and fisco = 9.0(±0.6) x 10"^. Figure M shows the locations of the 
CO emission in gray shade on top of the temperature and density profiles (top and middle 
panels) and the two vertical boundaries {ag = 0.79 and Tco = 19K) at R = 200 AU (bottom 
panel) . 

Figure [11] shows the effect of the CO freeze-out temperatures Tco on the radial column 
density of ^'^CO in the model. Without any freeze-out, the ^^CO column density follows the 
exponential taper of the H2 density profile. The freeze-out at 19 K leads to a significant drop 
in the gas-phase ^^CO column density beyond a radius of ~155 AU (or 310 AU in diameter), 
which we directly resolve. Even though the ^^CO column densities are different by orders 
of magnitude , the line emission difference projected to the line-of-sight from the outer disk 
is very subtle and high spectral resolution data is fundamental for resolving any molecular 
abundance structure changes at radii beyond ~155 AU. Indeed a similar analysis of the C^^O 
J=2-l emission alone indicates that the data can be fit equally well (or perhaps better) by 
models that do not include CO freeze-out (i.e. the lower boundary is the midplane, z = 0). 
However, this apparent inconsistency is likely the result of the 10 x lower velocity resolution 
of the C^^O data. Therefore, we don't fit for CO vertical boundaries from the emissions of 
C^^O 2-1 and C^^O 3-2 due to their weaker signals and limited spectral resolutions. 

Keeping that resolution effect in mind, we adopt the abundance boundaries derived 
from the ^"'CO 2-1 data and then fit the CO, C^^O, and C^^O data with only the fractional 
abundances as free parameters. For the as and Tco derived above, we find the fractional 
abundances of CO, C^^O and C^^O to be 6.0(±0.3) x 10-^ 1.35(±0.20) x 10"^ and 3.5(±1.1) x 
10-s respectively, corresponding to CO/^^CO = 67±8, CO/C^^O = 444 ±88 and CO/C^^O 
= 3.8 ± 1.7. Our de rived isotopic ratios are all consistent with the quiescent interstellar gas- 



phase values, which IWilsonI (119991 ) find in the local ISM to be CO/^^cQ = 69±6, CO/C^^O 
= 557 ± 30, and CO/C^^O = 3.6 ± 0.2. Our final model parameters are listed in Table HI 
The models are directly compared with the data channel maps in Figures [12] and Figure [12] 
(with the velocities binned in 1 km s~^ channels). Table H] also shows the best-fit fractional 
abundances in a model of HD 163296 that does not include CO freeze-out; the ^^CO/C^^O 
ratio is determined to be 30.6 ± 6.1, about four times higher than the value 6.7 ± 1.4 derived 
from our best-fi t models tha t consider CO freeze-out at 19 K (Table Hj) or the local ISM 



value 8.1 ± 1.1 flWilsonlll999l ). This provides an indirect evidence that we should also take 
into account CO freeze-out in the C^'^0 2-1 data analysis. 
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5. Discussion 

5.1. Modeling Dust Emission and CO Line Emission 

We have modeled multiple emission lines of CO and its isotopologues from the disk 
around HD 163296 in the context of an accretion disk model structure, grounded in obser- 
vations of the broadband SED and resolved millimeter continuum emission. The goal of 
this modeling effort is to develop a more consistent examination of the connection between 
the gas and dust phases in the disk. While our modeling framework does not treat the 
complexity of a completely self-consistent, simultaneous description of the energy balance 
and chemistry between the gas and dust, it effectively employs a set of parameters that can 
retrieve molecular abundance information in a way that captures the essential character of 
the layered disk structures predicted by those more sophisticated models. Most importantly, 
our approach directly addresses two common issues that are noted in much simpler structure 
models: (1) a radius (or size) discrepancy between the dust and CO emission, and (2) the 
degeneracy in the vertical temperature structure for models based solely on the dust emission 
(i.e., the SED). 

A longstanding problem with disk models has been the seemingl y different radial dis - 



tributions of dust and gas when each is considered independently ( JDutrey et al.l 120071 ). 



Isella et al.l ( 120071 ) presented multi-wavelength millimeter continuum and CO isotopologue 
observations of the disk around HD 163296 and found a significant discrepancy between 
the outer radius derived for the dust continuum (200 ± 15 AU) and t hat derived from CO 



emission (540 ± 40 AU) in their truncated power law models. However, iHughes et al.l (120081 ) 
showed that models with tapered outer edges can naturally reconcile the apparent size dis- 
crepancy in dust and gas millimeter imaging. The successful fitting of CO isotopologues in 
the HD 163296 disk using the SED-based disk model with an exponentially tapered outer 
edge, without invoking an unknown or unconstrained chemical effect, provides new support 
for the necessity of including this feature in outer disk structure. 

The SED alone does not provide any direct information on the temperature structure 
of the intermediate disk layers where the CO (and other molecular) emission is generated. 
In our modeling framework for the HD 163296 disk, the SED (and millimeter continuum 
images) can be fit equally well for a wide range of vertical temperature/density profiles, 
highlighting the degeneracy of the dust data with the parameter Zbig, the height marking the 
transition between the small grains in the disk atmosphere and the big grains concentrated 
toward the midplane. However, we have found that resolved observations of optically thick 
CO lines at a range of excitations can be used to place stringent constraints on the vertical 
temperature structure. Previous analysis of the CO J=2-l and 3-2 lines from the HD 163296 
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disk suggested that gas temperatures were always higher than 20 K , ruhng out CO freeze-out 
as a cause for the depletion of CO abundances (llsella et al.l 120071 ). But as we discussed in 
§4.2 (see Figure [H]), the small excitation leverage between those low- lying transitions is not 
a strong discriminant of the temperature profile. Here, we make use of the higher-excitation 
J=6-5 line to better constrain the vertical structure of the disk, and find a colder midplane 
that is consistent with significant CO freeze-out. Although observations of these various CO 
transitions is expensive, there are likely other molecules that emit at nearby frequencies and 
can be used to trace a sufficient range of excitation conditions (e.g., HCO^). The essential 
point is that the temperature structure in the intermediate layers of a disk can only be 
constrained well using several optically thick emission lines that probe a range of excitation. 



5.2. The "CO Snow Line" 



Theoretical models of the molecular abu ndances in protoplanetary disks predict three 
distinct vertical layers (see lBergin et al.ll2007l ). At large heights in the disk atmosphere, tem- 
peratures are high and molecular abundances are comparable to those in a standard photon- 
dominated region (PDR). At intermediate heights, between the midplane and atmosphere, 
warm temperatures are suitable for high gas-phase abundances of many molecules. And at 
the deepest layers near the midplane, temperatures are low enough that many molecules are 
deplet ed from the gas phase and frozen onto the mantles of dust grains. lAikawa &: Nomura 
(120061 ) and others, have used this kind of layered physical structure in their chemical reac- 
tion network models and shown that CO is expected to be abundant in the warm layer, with 
/co ~ 10^"^ almost independent of radius. The boundaries of that layer are defined by the CO 
freeze-out in the midplane and the photodissociation of CO by high-energy stellar photons in 
the disk atmosphere. The modeling of the SMA observations of the HD 163296 disk provides 
direct support for this basic model structure. Our estimate of the upper boundary as is in 
excellent agreement with models of the upper PDR layer, where CO is photodissociated into 
its ator aic constituents by stellar UV and/or X-ray radiat ion at column densities < 10^^ 
cm~^ |3| (JAikawa fc Nomurall2006l : iGorti fc Hollenbachll2008l ). The lower boundary that best 
fits the ^^CO data, corresponding to a temperature of 19 K, is also in excell ent agreement 



with the laboratory studie s of CO freeze-out temperature onto dust gr ains (ICoUings et al. 



Hersant et al. 



p et aJ 
l2009h 



20031: iBisschop et al.M2006l ). The main effect of CO photodesorption (lOberg et al. 



2007 



and other related non-thermal ice desorption mechanisms is to push the 



■^A total gas column of 10^^ cm~^ corresponds to 1 mag of extinction (A^, =1) at visible wavelengths 
when the gas-to-dust ratio is 100. But in our model the CO photodissociation front is located at Ay ^1 
because of the large dust depletion factor as measured with e = 0.003 in the upper disk layers. 
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lower boundary to colder temperatures. Since our best-fit CO lower boundary is consistent 
with thermal desorption, there is no evidence for efficient CO photodesorption in the disk of 
HD 163296. 

Since the term "snow line" is often discussed in planetary formation as some point where 
the temperature in the midplane would drop below the water ice sublimation level, e.g. in 
the "minimum-mass solar nebular" model prescribed by iHayashil (119811 ) , we propose to use 
the term "CO snow line" to indicate where CO freezes out in the disk. It was generally 
believed that the dust temperature in the disks of Herbig Ae stars is high enough that even 
the temperature close to the disk midpla ne is above the temperature of freeze-out of CO 
(JDutrey et al.l 120071 : iJonkheid et al.l 120071 ) . In this work, the primary line of evidence for 
CO freeze-out in the disk of HD 163296 comes from the detection of sharp reduction of the 
effective ^'^CO J=2-l emission area in the outer disk. This is demonstrated in Figure [TUJ 
where the ^^CO visibility is best fit with the lower boundary at 19 K, indicating no or very 
few ^^CO emission from the midplane area below this boundary, as shown in Figure [H The 
CO snow line marks drop off at least two orders of magnitude in the column density of ^^CO 
beyond 155 AU, which is clearly resolved by our interferometric observations. However, the 
emission lines from CO and its isotopologues are still detectable out to large radii, e.g. >500 
AU, which indicates that there are still substantial amount of CO toward surface layer of 
the outer disk where the temperature is higher than 19 K. 

The presence of this CO snow line has a direct impact on some key elements of chem- 
ical networks in a disk, especially related to the ionization fraction. The CO abundance 
is expected to be linearly correlated with the HCO+ abundance, but has a strong anti- 
corre lation with N2H"^ sinc e the main destruction pathway for the latter is in reactions with 
CO (jj0rgensen et al.ll2004l ). In the midplane, the depletion of CO beyond the snow line 
should enhance the deuterated molecular ions, H^D"*", D^H"*" and Dt, such that they become 
the most abundant ions in the midplane ( ICeccarelli fc DominikI 120051 ). So, resolving the CO 
snow line will eventually help constrain the ionization fraction and the extent of deuterium 
fractionation in the disk interior. 



5.3. CO Isotope Ratios 

The model developed here indicates that the gas in the intermediate layer of the 
HD 163296 disk has CO isotopologue abundance ratios that are consistent to those found 
in the molecular interstellar medium. That finding is in contrast with some previous stud- 
ies of this and other disks th at instead sugges t ed an increasing amount of depletion from 
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CO to i=^CO to C^^O (e.g., Ilsella et al.l l2007l : iDutrev et al.l Il994 |l996|). That difference 
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is likely a manifestation of the underlying modeling approach, as well as the constraining 
power now available to us from the high-excitation J=6-5 line and the rare isotopologues. 
When we include CO freeze-out, our models also show no evidence for radial variations of 
the CO abundance with respect to the H2 densities inferred (indirectly) from the dust emis- 
sion. Previous studies have found (much) steeper CO surface density profiles compared to 
dust (e.g., Pietu et al. 2007, Dutrey et al. 2008). To compare with those findings, we have 
used an alternative power- law prescription for the radial abundance fractions, fi{r) oc r^*, 
and attempted to re-fit our data. The resulting best-fit values of Pi for ^^CO and C^^O are 
0.0 ± 0.1 and 0.1 ± 0.2, respectively: both consistent with a radially constant abundance 
profile. So, when taking into account CO freeze-out at low temperatures, there is no clear 
evidence that the radial distribution of CO deviates from the underlying total gas surface 
density. The steep slopes inferred in previous work may be a manifestation of the sharp 
abundance drop beyond the CO snow line. 

This model demonstrates that, given sufficiently strong constraints on the vertical tem- 
perature gradient, a unified model of the gas and dust disk can be constructed that conforms 
to standard ISM-based assumptions about the isotopic abundances and the CO chemistry. 



Visser et al.l ( 120091 ) present a photodissociation model for CO isotopologues including 
newly updated depth-dependent and isotope-selective photodissociation rates. They find 
grain growth in circumstellar disks can enhance iV(^^CO)/A^(C^'''0) and A^(^^CO)/A^(C^®0) 
ratios by a factor of ten relative to the initial isotopic abundances. We have fit the cr^ and 
fractional abundances for C^^O and C^®0 and found that different pairs of ag and f(C^^O) or 
f(C^^O) can fit the data equally well, i.e. we can not constrain ag due to the limited signal-to- 
noise ratios of C^^O J=2-l and C^^O 3-2 data. Future observations with greater sensitivity 
in the CO isotopologues will be essential to investigate the isotope-selective photodissociation 
on the distributions of CO isotopologues and hence the local isotopic ratios at different layers 
of the disks. 



5.4. Comparison of HD 163296 and TW Hya 



In a series of papers, lOi et al.l ( 12004 |2006| ) modeled the CO J=2-l, 3-2, and 6-5 line 



emission from the disk around the cooler young star TW Hya. There, the best fit model to 
the optically thick CO J=2-l and 3-2 lines tended to underestimate the CO 6-5 emission. 
To fit all three CO lines simultaneously required additional heating of the disk surface, 
suggested to be the result of the intense X-ray irradiation field for that source. No such 
additional heating is required to explain the HD 163296 data. On the contrary, the CO 
J=6-5 emission required a lower disk interior temperature compared to those predicted by 
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the typical accretion disk models based on the SED, forcing us to modify the Zbig parameter 
to fit the data. 

There are at least two plausible reasons for the different vertical structures we infer 
for these disks. First, is the high-energy radiation field from the two stars. Chandra X-ray 
observations reveal a point-like object within Cf.^25 (30 AU) of HD 163296, with an X-ray 



lumi nosity of 4 x 10 er gs s (jSwartz et al.ll2005l ). about 5x lower than that from TW 



Hya (JKastner et al.lll999l ). Obviously the suggestion that X-ray heating is dependent on the 
stellar type will require a much larger sample to confirm, as well as more detailed radiative 
transfer modeling of the high-energy heating of these disks. However, it is an exciting 
prospect that the interior structures of these disks may be very different depending on the 
connection between the stellar mass and the irradiation environment, with implications both 
for disk evolution and chemistry. A second possibility may lie with the level of turbulence 
in these disks. In our model, the vertical structure is governed by the settling of large dust 
grains. The timescale for gravitational sedi mentation is expected to be sho rt, with big grains 
reaching the midplane in less than ~1 Myr ( iDuUemond fc Dominikll2004j ) unless the dust is 
effectively stirred up by turbulence. The more vertically extended population of big grains 
required to model the HD 1 63296 disk might be commensurate with more efficient stirring, 
or rather inhibited settling. iHughes et al.l ( 120111 ) suggest a 300 m s~^ turbulent linewidth 
in the CO layer of the HD 163296 disk, much larger than the <40 m s~^ they find for the 
TW Hya disk. Perhaps turbulent stirring is responsible for lofting the big grains above the 
midplane in the HD 163296 disk, enhancing the radiative cooling. 



6. Summary 

We present CO line and dust emission observations of the disk around the Herbig Ae star 
HD 163296 and develop a model framework that describes in a consistent way the spectral 
energy distribution, resolved millimeter dust continuum data, and multiple emission lines 
of CO and its isotopologues. The fitting results support the general picture of CO vertical 
distribution regulated by photodissociation at a surface where the hydrogen column density 
is < 10^^ cm~^ and by CO freeze-out at depths below 20 K in the midplane. The main 
conclusions are summarized here: 



1. We confirm a tapered exponential edge in the surface density distribution, an outcome 
of the similarity solution of the time evolution of acc retion disks, can ac count for the 
size discrepancy in dust continuum and CO emission (JHughes et al.ll2008[ ). This means 
that constraints from dust modeling can be applied to CO gas modeling, significantly 
reducing the number of free parameters. 
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2. We find in the disk model the transition between the "small" and "big" dust grain pop- 
ulations, z = zug in units of the gas scale height H, regulates the vertical temperature 
and density profiles between the disk midplane and surface. Multiple transitions of 
CO, especially the CO 6-5 line, which requires much higher excitation, can be used to 
constrain the location of zug, thus the temperature structure in the disk intermediate 
layer. We find the resulting disk model for HD 163296 has a cold midplane populated 
by large grains with a large scale height (^big = 2H). Since the vertical temperature 
in the HD 163296 disk is governed mainly by the settling, or lack of settling, of large 
grains, a possible explanation for th e relatively cold inte rior of HD 163296 disk may 



be a high level of gas turbulence (see iHughes et al.ll201ll ) 



3. Using the model with temperature structure constrained by the CO multi-transition 
analysis, we fit for the location of emission from ^'^CO and constrain the vertical distri- 
bution of the emission region to lie between a lower boundary set by the temperature 
where CO freezes out onto dust grains (at heights where T < 19 K) and an upper 
boundary where densities are low enough that stellar and interstellar radiation can 
photodissociate the CO molecule (where the vertically integrated hydrogen column 
density from the disk surface is < 10^^ cm~^). The CO freeze-out produces a signif- 
icant drop in the gas-phase CO column density beyond a radius of about 155 AU, 
effectively a CO snow line that is resolved directly by the observations. The CO deple- 
tion, generally found in disks, can be successfully accounted for considering both the 
CO freeze-out and photodissociation. 

4. Taking the CO freeze-out into consideration in the disk model and assuming CO iso- 
topologues sharing the same spatial distribution, the isotopic ratios of ^^C/^^C, ^^O/^^O 
and ^^O/^^O are consistent with the standard quiescent interstellar gas-phase values 
and show no evidence for unusual fractionation. More sensitive data is needed to inves- 
tigate the distribution differences among those CO isotopologues and the local isotopic 
ratios at different layers of the disks. 

This detailed investigation of the HD 163296 disk demonstrates the potential of a staged, 
parametric technique for constructing unified gas and dust structure models and constrain- 
ing the distribution of molecular abundances using resolved multi-transition, multi-isotope 
observations. The analysis provides the essential framework for more general observational 
studies of molecular line emission in protoplanetary disks. 
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Table 1. Observational parameters for SMA HD 163296 observations 



Parameters CO 2-1 CO 3-2*^ CO 6-5 ^^CO 2-1 C^^O 2-1 C^^O 3-2 



Rest Frequency (GHz): 


230.53800 


345.79599 


691.47308 


220.39868 


219.56036 


337.06113 


Observations'' 


2010 May 17 (C) 


2009 May 6 (C) 


2007 Mar 20 (C) 


2010 May 15 (C) 


2010 May 15 (C) 


2005 Aug 23 (C) 




2010 Sep 14 (E) 


2009 Aug 22(E) 




2010 Sep 11 (E) 


2010 Sep 11 (E) 




Tsys 


93-207 


147-273 


1368-4764 


82-146 


82-146 


210-603 


Baselines 


6-151 


10-146 


23-158 


6-145 


6-145 


8-80 


Spectral resolution (km s"'^): 


0.066 


0.044 


0.88 


0.069 


0.55 


0.72 



''C: compact configuration; E: extended configuration. 
''Observational parameters from iHughes et al.l pOllh 



to 

to 
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Table 2. Continuum and emission line results 



A(mm) 


Beam 


PA 


Flux (Jy) 


1.36 


3'.'5 X 2'.'2 


55.8° 


0.615±0.004 


1.33 


4'.'7 X 2'.'8 


24.2° 


0.670±0.007 


1.11 


3'.'0 X I'.'g 


57.8° 


1.03g±0.007 


0.88 


1'.'8 X 1'.'3 


47.6° 


1.74±0.12 


0.44 


2'.7 X 2'.'2 


-6.g° 


7.5±0.5 


Transitions 


Beam 


PA 


Integrated Intensity^ (Jy km"^) 


CO 2 1 


2'.'1 X 1'.'8 


10.5° 


54.17±0.3g 


CO 3-2 


1'.'7 X 1'.'3 


46.7° 


g8.72±1.6g 


CO 6-5 


2'.7 X 1'.'4 


1.8° 


58.66±6.44 


i^CO 2-1 


I'.'g X 1'.'8 


5.5° 


18.76±0.24 


C^^O 2-1 


I'.'g X 1'.'8 


5.5° 


6.30±0.16 


C^^O 3-2 


3'.'2 X 2'.'2 


13.1° 


11.64±0.76 



^Intensity averaged over the whole emission area. 
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Table 3. Physical model for the disk of HD 163296 



Parameters 



Values 



Stellar and accretion properties 



Spectral type 

Effective temperature: T^.(K) 
Visual extinction: A^ 
Estimated distance: d(pc) 
Stellar radius: R*(R0) 
Stellar mass: M,,(Mq) 
Accretion rate: M{Mq yr^^) 



Al 

9333 

0.3 

122 

2 

2.3 
7.6x10"^ 



Disk structure properties 



Disk mass: Md(M0) 

Characteristic radius: Rc(AU) 

Viscosity coefficient: ao 

Depletion factor of the atmospheric small grains: e'^ 

Inner wall radius: R^aii(AU) 
Inner wall scale height: z^„a//(AU) 



0.089 
150 
0.019 
0.003 
2.0 
0.6 
0.1 



Disk geometric and kinematic properties 



Inclination: i(deg) 
Systemic velocity: VLS'/j(km s~^) 
Turbulent line width*^: Vtu^b(km s^^^ 
Position angle: P.A.(deg) 



44 ± 2 
5.8 ± 0.2 

0.2 
133 ± 3 



'^See definition in paper or lD'Alessio et al.l ( l2006l ). 
''Gas scale height. 
'^Fixed parameter. 



Table 4. Fitting results: fractional abundances and distributions 



Parameters CO 2-1 & 3-2 ^^qq 2-1 C^^O 2-1 C^^O 3-2 

Midplane freeze-out temperature (K) 19.0* 19.0 ± 0.3 19.0* 19.0* 

as 0.79* 0.79 ± 0.03 0.79* 0.79* 

Fractional abundance (6.0±0.3) xlO^^ (9.0±0.6) xlO"^ (1.35±0.20) xlO"^ (3.5±1.1) xlO^^ 

Fractional abundance (no freeze-out) (6.5±0.3) xlO'^ (4.6±0.3) xlO^^ (1.5±0.2) xlO^^ (7.0±2.2) xlO^^ 

*Parameter values adopted from ^'^CO 2-1 fitting. 



to 
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l""l""U"'l""l""l 



■10 -5 5 10 15 20 

Vlsr (km s"') 



^CO 2-1 




l""l""U"'l""l""l' 



■10 -5 5 10 15 20 

Vlsr (km S"') 



Fig. 1. — Spatially integrated spectra at the peak continuum (stellar) position of HD 163296. 
The fluxes are averaged over the emission areas (4" box for CO 6-5; 10" box for other lines). 
The vertical dotted line indicates the position of the fitted Vlsr- 
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Fig. 2. — Channel maps of the CO J=2-l hne emission from the disk around HD 163296. 
The LSR velocity is indicated in the upper right of each channel, while the synthesized beam 
size and orientation (2'.'1 x I'.'S at a position angle of 10.5°) is indicated in the upper left 
panel. The contours are 0.18 Jy Beam~^ (la) x [3, 6, 9, 12, 15, 18, 21, 24, 27]. The star symbol 
indicates the disk center. 
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Fig. 3. — Same as Figure [21 but for ^"^CO 2-1. The beam size is 1'.'9 x I'.'S at a position 
angle of 5.5° and the contours are 0.08 Jy Beam"! (la) x [3, 6, 9, 12, 15, 18, 21, 24]. 
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Fig. 4.— Same as Figure El but for C^^O 2-1, C^^O 3-2 and CO 6-5. For C^^^O 2-1, the 
beam size is 1'.'9 x I'.'S at a position angle of 5.5° and the contours are 0.025 Jy Beani"^ (la) 
X [3, 6, 9, 12, 15, 18, 21, 24]; for C^^O 3-2, the beam size is 3'.'1 x 2'.'3 at a position angle of 12.8° 
and the contours are 0.2 Jy Beam"^ (Icr) x [2, 3, 4, 5]; for CO 6-5, the beam size is 2'.'7 x 2'.'2 
at a position angle of -5.6° and the contours are 3.0 Jy Beam^^ (la) x [2, 3, 4, 5, 6, 7]. 
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Fig. 5. — Top panels (Left): The predicted continuum image vs data at 271 GHz. Contours 
are shown at Sa intervals. (Right): The predicted 271 GHz continuum visibihty profiles 
for Models with Rc=l25 AU(dotted line), 150 AU(solid line), 175 AU(dashed line). Bottom 
panels: The predicted continuum image and visibilities vs data at 341 GHz. 
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Fig. 6.— T he H P 163296 SEP . Optical and in fr ared bands: U BVJHKLM from 



Malfait et all (119981 ): U BVRl from iHillenbrand et all (|l992|); BV from lOudmaiier et al. 



J200lh: BVR I JHK from IXannirkulam et al.l J2008I) 



Berrilli et al. (11992): JHK from 



spect ra (jAcke fc van den 



Ancker 



^iroa et al. (200 



2004 



Thi et al. 



JHKLM and 8-13 /imm bands from 
) an d 2MASS JCutri et al.l [2OO3I ): ISO 
200 ih: the CTIO and Keck telescopes 



data ( iJayawardhana et al.l 1200 ll ): IRAS (JBeichman et al.lll988[ ): Spitzer-IRS (retrieved from 



Spitzer archive we bsite: irsa.ipac.calt ech.edu). Milhmeter fluxes are from this paper (marked 
with circles) and llsella et al.l ( l2007l ). Overlaid on the SEP are the models which have: 
M = 7.6 X 10-s Mq/vts, ao = 0.019, e = 0.003, R^ = 150 AU, cosi = 0.72, a vertical wall at 
the dust sublimation radius, calculated assuming Tg^b = 1500 K, with Rwaii = 0.6 AU and 
a height z^aii = 0-^7Rwaii- The stellar photosphere (the shor t -long -dashed line) is from the 
library for population synthesis from iBruzual A. fc Chariot! ( 119931 ). The model SEP lines 
are laid out for models with different values of z^g/H = 0.5(long-dashed), l.O(short-dashed), 
1.5 (dot-dashed), 2.0(sohd) and 2.5 (dotted). 
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Fig. 7. — Upper panel: Temperature vs height at i? = 200 AU for disk models with M = 
7.6 X 10^^ M^/yrs, Uq = 0.019, e = 0.003, R^ = 150 AU, and different values of z^g/H = 
0.5, 1, 1.5, 2.0, 2.5, 3 (from left to right). The fiducial model, with zug = 2iJ is showed with 
a solid line. Lower panel: Density vs height at i? = 200 AU for the same models. The lines 
are laid out from top to bottom for models with different values of zug/H = 0.5, 1, 1.5, 2.0, 
2.5, Sat z> 100 AU. 
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Fig. 8. — CO J=2-l, 3-2, and 6-5 model spectra (black lines) overlaid with the HD 163296 
spectra in grey shade. The left, middle and right column are the simulated models with 
zugl H = 1.5, 2.0, 2.5. Note that all the models have been sampled at the same spatial 
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Fig. 9. — The irradiated accretion disk structure model for HD 163296 with Z\^\JH = 2.0 
and i?c= 150 AU. Top and middle panels: Temperature and density profiles are indicated 
by solid and dashed lines, respectively. The CO emission area constrained by the best-fit 
vertical boundaries is shown in grey shade. Bottom panel: The vertical distributions of 
temperature and cr^ at R = 200 AU are indicated by solid and dashed fines, respectively. 
The vertical dotted lines show the best-fit locations of the lower boundary (Tco = 19 K) 
and upper boundary (as=0.79) for CO. 
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Fig. 10. — Iso-x^ surfaces of CO freeze-out temperature (Tco), CO fractional abundance 
(fco) and as for ^^CO. Contours correspond to the 1-5 a errors. 
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Fig. 11. — Models of ^^CO radial column densities for no freeze-out (dotted), freeze-out at 
18 K (dashed) 19 K (solid) and 20 K (dot-dashed). Each model has been scaled to fit the 
^^CO 2-1 emission. Note that the CO snow line is at a radius of 155 AU for Tco = 19 K 
and it increases from 135 to 175 AU when Tco decreases from 20 to 18 K. 
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Fig. 13. — For each panel, the top rows are the velocity channel maps of the CO 2-1, 
3-2 and 6-5 emissions toward HD 163296, respectively (velocities binned in 1 km s^^). 
Contours are 0.15 Jy Beam-^ (Ict) x [3, 6, 9, 12, 15, 18, 21, 24] for CO 2-1; 0.2 Jy Beam"! 
(la) X [3, 6, 9, 12, 15, 18, 21, 24, 27] for CO 3-2; 3.0 Jy Beam"! (la) x [2, 3, 4, 5, 6, 7, 8] for CO 
6-5. The middle rows are the best-fit models and the bottom rows are the difference between 
the best-fit models and data on the same contour scale. 



